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In the present paper, we investigate the dynamics of a Bose-Einstein condensates (BEC) loaded 
into an deep optical lattice of ID, 2D and 3D, both analytically and numerically. We focus on 
the self-trapping state and the effect of the system dimension. Under the tight-binding approxi- 
mation we obtain an analytical criterion for the self-trapping state of BEC using time-dependent 
variational method. The phase diagram for self-trapping, soliton, breather, or diffusion of the BEC 
cloud is obtained accordingly and verified by directly solving the discrete Gross-Pitaevskii equation 
(CPE) numerically. In particular, we find that the criterion and the phase diagrams are modified 
dramatically by the dimension of the lattices. 
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I. INTRODUCTION 



Recently, the Bose-Einstein condensates (BECs) in an 
optical lattice attract much attention both experimen- 
tally and theoretically [1]. This is partly because the 
physical parameters, e.g., lattice parameters, interaction 
strength, etc., in such a clean system ("clean" means per- 
fect periodicity without impurity) can be manipulated 
at will using modern experimental technique like Fesh- 
bach resonance. Taking this advantage people have ob- 
served well-known and long-predicted phenomena, such 
as Bloch oscillations [2] and the quantum phase transition 
between superfluid and Mott insulator [3]. More impor- 
tantly, there are new phenomena that have been either 
observed or predicted in this system, for example, non- 
linear Landau-Zener tunneling between Bloch bands [4, 5] 
and the strongly inhibited transport of one-dimensional 
BEC in an optical lattice [6]. 

Another intriguing phenomenon, namely, self- 
trapping, was recently observed experimentally in this 
system [7], In such a experiment, it was observed that a 
BEC cloud with repulsive interaction initially loaded in 
an optical lattice was self-trapped or localized and the 
diffusion was completely blocked when the number of 
the degenerate bosonic atom was larger than a critical 
value. This is quite counterintuitive. Even without 
interaction, a wave packet with a narrow distribution 
in the Brillouin zone expands continuously inside a 
periodic potential. When there is an interaction between 
atoms and it is repulsive, one would certainly expect 
the wave packet to expand. On the other aspect, the 
self-trapping phenomenon also attract many theoret- 
ical efforts and a lot of meaningful results have been 
obtained[8-14]. However, most of these analyses are 
limited to ID system. Beyond these ID approximation, 
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the localized states[15]-[16] , self-trapping and other 
quantum phenomena[17]-[20] have been observed in 
high-dimensional optical lattices. To our knowledge, 
however, there are no systematical theoretical analyses 
on the dynamics of BEC loaded into multi-dimensional 
optical lattices. There leaves some important problems, 
e.g., what is the effect of the system dimension on the 
dynamics of BEC in optical lattices ? 

In this paper, we address this issue both analyti- 
cally and numerically focusing on the intriguing self- 
trapping phenomenon, within the tight-binding approx- 
imation. Based on the discrete GPE and exploiting a 
time-dependent variational approach, we obtain an an- 
alytical criterion for the self-trapping state of BEC and 
the phase diagram for self-trapping, soliton, breather, or 
diffusion of the BEC cloud. Our results show that, the 
change on the dimension of the system leads to many 
interesting consequence: i) Self-trapping state can also 
emerge in higher-dimension systems, however, compared 
to the ID situation its parameter domain in the phase 
diagram is largely shrunken; ii) The stable moving soli- 
ton and breather solutions of the wavepackets that exist 
in ID system, no longer stand in 2D and 3D systems, 
iii) In addition, when the self-trapping occurs, the ra- 
tio between the final BEC cloud width and the initial 
wavepacket width is greatly modified in 2D and 3D sys- 
tems compared to the ID case. 

Our paper is organized as follows. In Sec. II, under the 
tight-binding approximation, we obtain a discrete GPE 
that governs the dynamical behavior of BEC in deep 
31? periodic optical lattices. In Sec. Ill, using the time- 
dependent variation of a Gaussian profile wavepacket, the 
dynamics of BEC is discussed analytically and the crit- 
ical conditions for occurrence of self-trapping state are 
obtained. Moreover, the phase diagram for self-trapping, 
soliton, breather, or diffusion of the BEC matter waves is 
obtained accordingly. In Sec. IV, directly solving GPE is 
exploited to confirm our theoretical predictions. Sec.V, 
are our discussion and conclusion. 
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II. NONLINEAR DISCRETE MODEL 

We focus our attention on the situation that a BEC 
is loaded into a 3D optical lattices. In the mean field 
approximation, the 3D BEC dynamics at T=0 satisfies 
the GPE 

where VP is the wave function of the condensate. The 
coupling constant Airti 2 a s /M, which relate to the mass 
of the atoms M and the repulsive (attraction) s-wave 
scattering length a s > 0(a s < 0), is the nonlinear 
coefficient corresponding to the two body interactions, 
Vopt(~?) — U max (cos(2Kx)+ cos(2kj/) + cos(2k2;)) is the 
external optical potential created by three orthogonal 
pairs of counter-propagating laser beams [21], k is the 
wave number of the laser lights that generate the optical 
lattice. All the variables can be re-scaled to be dimen- 

sionlcss by the system's basic parameters t ~ tv 1 
* ~ -7=*, then GPE (1) becomes 

»flr = - va * +v ( r )* + ^-i*i a * w 

where N is the total number of atoms, Vi^r) = 
Vb[cos(2x)+ cos(2y) + cos(2z)], the parameter Vb = 
V max /2E]i is the strength of the optical lattices, with 
Eji = h 2 K 2 /2M is the recoil energy of the lattice, v = 
Er/H. The associated Hamiltonian is 

h= J a-? + v{-?m 2 + 2 ^^V i 4 ] 

(3) 

In order to understand the dynamical behaviors of 
BEC in deep 3D optical lattices, it is worthwhile to study 
the system in certain extreme limits. The particular limit 
that we now insight into is the tight-binding approxima- 
tion, where the lattice is so strong that the BEC system 
can be looked as a chain of trapped BECs that are weakly 
linked [8]. Under the tight binding limit, the condensate 
parameter can be written as 

*(~r,*)= ^n,m,k{t)<i>^r - "r„ (mife ) (4) 

n,m,k 

where Vn,m,fe(*) = \/N 1hm ^(t)e ien - m - k(t) is the 
(n, m, fc)th amplitude, N njTn ^ and 6 n _ m ^{t) are the 
number of particles and phases in the well (n,m,k). 
J2 n ,m,k l*n,m,fc| 2 = 1 , 4>(j - ^n,m,k) is the conden- 
sate wavefunction localized in the well (n, to, k), with the 
normalization J d'r'(j) 2 lmk = 1, and the orthogonality 
of the condensate wave functions J d~r<t>n,m,k<t>n+i,m.,k — 

0,J d~T '4> n , m ,k<Pn,m+l,k — 0, / d~t 4>n,m,k4>n,m,k+l — 0. 

Inserting the nonlinear tight-binding approximation (4) 



into the GPE (2) and integrating out the spatial degrees 
of freedom, following re-scale the time as t — > 2Kt, we 
find the reduced GPE 



i g 7 "' = --^(lpn-l,m.k + i>n+l,m,k + i>n,m-l,k 

+1pn,m+l,k + V'n.m.fe-l + V'n.m.fc+l) 

+ (A + A\tpn,m,k\ 2 )^n,m,k (5) 

where A = ^ / [(V0„, m , fc ) 2 + V{~r )0„, m , fc 2 ], 

A = A ^Mfi^ Jd?4>n,m,k 4 is the ratio of the 
nonlinear coefficient, induced by the two-body in- 
teratomic interactions. And the included parameter 

K — ~ J dr{(V<j}n+l, m ,k^4>n, m ,k + V0 nim+ l,fcV0 n , mi k + 

V0„, mi fc + iV0„, m ,fc) + V r [cos(2a;)^„ + i ;mi fe(/>„ iTO;fe + 
cos(2y)(/)„ !m+ i ifc 0„^ m!fc + cos(22;)0 n;mife+ i0„ iTO;fe ]} is the 
tunneling rates between the adjacent sites (n, to, fc) and 
(n+l,m+l,fc + l). In order to make an agreement 
with the experiments done with s7 Rb, the two-body 
interaction term is expected to be positive, that is A > 0. 

Equation (5) can be written as ip n m k — 
dH/d(itp* m k ), where H is the Hamiltonian func- 
tion 



n.m,k 

+V'n,m,feV'n,m+l,fe + V'n,m,kV'n,m+l,k 
+V ; n,m,feCm,fe + l + C, m ,feV ; n,m,fe+l) 

+^|^,m,fc| 2 + ||^n,m,fe| 4 ] (6) 



III. ANALYTIC RESULTS BASED ON THE 
TIME-DEPENDENT VARIATIONAL APPROACH 

To study the dynamics of the BEC wavepacket loaded 
into multi-dimensional optical lattices, we consider 
the evolution of a Gaussian profile wavepacket, which 

we parametrize as ip n ,m,k(t) — \/P ex ¥>{~ j^" 1 ~ 
(m j" )2 - + i\p x (n - £*) + p y {m - £„) + Pz(k - 

60] + %Mn - £.) 2 + *„(m - £,) 2 + <J,(* - C,) 2 }, 
define a = x,y or z, where £ a (t) and i? Q (i) are the 
center and the width of the wavepacket in the a 
direction, p a (t) and S a (t) are their associated mo- 
menta in the a direction, and p is a normalization 
factor. The wavepacket dynamical evolution can be 
obtained by a variational principle from the Lagrangian 

L = T,n,7n,k i (i ) n,m,kfn,m,k ~ V£, m ,fc V>n,m,fc) ~ H > with 

the equations of motion for the variational parameters 
Qi(t) = ta,R a ,Pa,5 a given by j- t j±- = where H 

is given by Eq.(6). In calculation, we can replace the 
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sums over (n, m, k) in the Lagrangian with integrals 
when R a is not too small [8]. In this limit, the nor- 
malization factor is obtained as p = 



W2 ^JJ Ra, 

from J^jTooJ-oo \lpn,m,k\ 2dndmdk = L AU 

these parameters lead to L = —V(R a ,£ a ) — 



2V^ J n Rp 

0=x,y,z 

where a a = 



cos p a 



i 



Risi 



and V{R a ,£J 
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of f f dndmdkAe 

' J — oo J — oo J — oo 

The equations of motion of collective variables £ a (t), 
R a (t), p a (i), S a (t) are given by 



trap, V(R, £) = 0, the momentum remains a con- 
stant p . Initially, we set £ = and S = 0. In this case, 
the effective Hamiltonian becomes 



H 



-D 



cospo • e a + 



A 



2^ D R D 



(11) 



Clearly, the dynamical properties of the system not 
only governed by the parameter A, but also influenced 
by the dimension D. As discussed above, the dynamics 
of the wavepacket are modified significantly by the sign 
of cospo- So two cases with cospo > and cosp < 
will be discussed respectively. 



Pa 



1 dV 

1 _„ 

= -e a smp a 

2 S 2 
= cosp a (^--^)e ' 

A 

+ v^ 3 ri n Rp 

f)=x,y,z 



Ra — 



2 8V 

R a dR a 



Ra 



S a e CTq cospo 



(7) 



A. cospo > 

Let us first consider the case with cosp > 0, following 
this case there exist two nonlinear phenomcnons. 

Self-trapping: Equation (9) shows that when the di- 
verging effective mass m* — ► oo, the group velocity of 
the wave packet v g — > 0, this corresponding to a self- 
trapping solution. Thus, for t — > oo one has R > Rm&x ? 

S — > oo and £ — > 0, H — > 2 ^d rD ■ Because the Hamil- 
tonian is conserved, H = Ho > 0, so when Ho > 0, we 
get the maximum value of the width 



For simplicity, we assume: R x — Ry — R z — 



Px 



Py 
1 

1W 



P, 5 X = Sy 



8, cr x 



R 2 S 2 



, £ x = £ y = £ z = £. Then, the effective 



Hamiltonian reduce to 



H = -De- a cosp + V{R,0 



A 



2^ D R D 



(8) 



where D = 1, 2, or 3 reflects the dimension of the system. 
It can be seen from Eq.(8), D plays an important role in 
the effective Hamiltonian. 

The group velocity of the wave packet is 



Rn 



= (- 



a/z? 



A 

Rn 



2Dy / 7T COSpo ' 



(12) 



Diffusion: In this case, for t — > oo, one has R — > oo, 

5^0 and £ — > i sinpo, as well as the effective Hamilto- 
nian H —> —D cospo- So the condition — cos < i?o < 
corresponding to the diffusive region of the system. The 
critical condition between self-trapping and diffusion is 



given by H 
one can get 



0, i.e., —Dcospo • e 2R o + 



= 0, 



dH 



tanp 



dp m 
with the inverse effective mass 



1 

to* 



o 2 h 



dp 



— = cosp 



• e 



(9) 



(10) 



It is important to note that, because this effective mass is 
related to the quasi-momentum, the system will exist rich 
dynamical phenomenons. Especially, when cosp < 0, the 
effective mass is negative and the system can exist the 
localized solution, i.e. soliton solution. For a untilted 



A c = 2y / 7r' D DRq cospo • e 



(13) 



The self-trapping phenomenon occurs at A > A c , and the 
diffusion occurs at A < A c . By observing the equation 
(12), one readily obtain 



Ro/Rmax. — (1 X^ 1 ^ 



(14) 



It is clear that when D=l, Eqs. (13) and (14) reduce to 
the corresponding results given in [8]. 
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cospo < 



When cospo < 0, the dynamics of the system (11) can 
be governed by (see Eq.(7)). 



-D| cos po| 



R, 



i? 4 



A 



^ D R D + 2 



R = ~7^ e CT | cos Po| 



(15) 



This system has moving soliton solution. The moving 
soliton solution corresponding to the fixed points R = 
and (5 = 0, £ =constant. One readily gets from (15): 




FIG. 1: The Phase diagram of Eq. (15) for 2D [(a), ft = 5, 
6o = 0] and 3D case [(b), ft = 10, 5o = 0] 



5 a =0,A C = 2V^ i? D_2 |cospo|e ™% (16) 

We now discuss the stability of this soliton state. Let 
5 = 5o + 5', R = Ro + R', and linearize Eq.(15) at the 
soliton state, we have: 



( S '\ _ A ( A{2-D-1/Rl)\(5'\ 

Kr'J 4 ^ D < +3 v- i? o )\R>) 

(17) 

one obtains its eigenvalue A 2 = \Rq{D - 2 + 1/Rq). It is 
clear that when D — 1, the soliton state is a center point 
(note that i? > 1), stable soliton solution and breather 
solution can exist in ID case. The dynamical behaviors 
of the wavepacket in ID case is discussed in detail in 
Ref. [8]. When D > 1, however, the soliton solution is a 
saddle point, no stable soliton and breather solutions can 
exist in the system. The numerical results for the phase 
diagram in the R — 5 space given by Eq.(15) for 2D and 
3D cases are shown in Fig.l. One can find that soliton 
and breather cannot exist in higher dimensional system. 
When A > A c , for t — > oo, one has R — > i? m i n < -Ro 
and 5 — -> oo. It means that the BEC wavepacket remains 
finite finally and the self-trapping phenomenon occurs. 
While A < A c , the situation is changed. It is shown 
that R — > oo and 5 — * 0, diffusion takes place. That is, 
the dynamics of the BEC wavepackets in ID, 2D and 
3D optical lattices are very different. In ID, soliton, 
breather, and self-trapping states can exist. In 2D and 
3D systems, however, only self-trapping state can exist 
when A > A c . 

In this case (cospo < 0), after the self-trapping occurs, 

R * Rmin 

< R , S — > oo and £ — > 0, H — > 2 r^n RD •> onc 

_ i 

can get from H — H = — D cosp • e + 

2-wtt rip. 



Ro/R min = (l + ^) 1/D (18) 

where A' = 2^/k D DR®\ cosp \e , and A > A c . 

The first and second panels in Fig. 2 show the analytical 
results for the dynamical diagram of the wavepacket in 
the cospo — A c plane. It reveals that the dynamics of the 
system can be influenced dramatically by the dimension. 
Compared with the ID case, the results are modified 
greatly in 2D and 3D systems. The region of the self- 
trapping decreases quickly with D, while the regions of 
the diffusion increases with D. We can clearly find that 
the soliton and breather solutions disappear in 2D and 
3D systems. It is also shown that, the effect of D on the 
regions of the self-trapping is more significant when the 
width of wavepacket is increased. 

The critical value A c against the initial wavepacket 
width Ro given by Eqs.(13) and (16) are also described 
in the third panel in Fig. 2. It is clear that the variation 
of self-trapping region against i?o is very different for ID, 
2D, and 3D cases, especially, for cosp < case. When 
cospo < 0, A c reaches a maximum in the beginning and 
then decreases with increasing Rq in ID. However, this 
relation is modified in 2D, where A c increases with Rq 
sharply when i? < 3 and finally A c remains an invari- 
able value when Ro > 3. Especially, A c increases with Ro 
linearly in 3D. Moreover, for a fixed Ro, A c in 3D case 
is several times, even one order of magnitude larger than 
that in ID and 2D. That is, if the initial wavepacket is 
given, the higher the dimension is, the greater the atomic 
interactions or atom number are needed in order to get 
into the self-trapping state. Thus it can be seen the di- 
mension has significant effect on the critical condition for 
the occurrence of the self-trapping state. 

Fig. 3(a) shows the ratio Ro/Rma* versus ratio A/A c for 
cospo > 0. When the self-trapping occurs, one can find 
that the final wavepacket width R max is larger than the 
initial wavepacket width Rq. For a fixed A/A c , Ro/R m nx 
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As is expressed, i? m i n is always smaller than Rq. More 
interestingly, i? m i n increases with Rq in ID and 3D, while 
-Rmin keeps a constant (0.75) and nearly independent with 
Rq. On the other hand, because the effect mass given by 
Eq.(10) is negative, the BEC cloud will finally squeeze 
into several lattice points for ID and 31? case and about 
two lattice points for 2D case. 




FIG. 2: The critical value A c against cospo (the first two 
panels) and initial width Ro (the third panel) 




FIG. 3: The value of Ro/R ma ,x for cospo > (a) and the 
minimum value -R m i n against the initial width Ro for cos po < 
0(b) 



increases with increasing D. It means that, when A > 
A c , the BEC cloud can easily get into self-trapping state 
with little adjustment of the initial wavepacket for higher 
dimension case. In addition, it is clearly that Ro/Rm^ 
increases with A/A c and close to 1 eventually for ID, 2D, 
and 3D cases. In other words, for sufficiently large A, the 
BEC cloud will quickly get into self-trapping state for all 
the case and the effect of dimension is weakened. 

Fig. 3(b) shows the minimum wavepacket width i? m i n 
against the initial wavepacket width Rq for cospo < 0. 



To confirm the above theoretical predictions, the di- 
rect numerical simulations of the full discrete GPE(5) 
are presented in Figs. 4-9. 

The numerical results for the dynamics of 2D 
wavepacket with small initial width (Rq = 2.5) are pro- 
vided in Fig. 4. It is explicit that when A > A c (first 
panel) the BEC cloud is quickly squeezed into two optical 
lattices finally and the width of the wavepacket remains 
finite (self-trapping). On the contrary, when A < A c 
(second panel), as it is expected the wavepacket will ex- 
panding indefinitely with time (diffusion). 

Figs. 5-7 show the numerical results of the wavepacket 
with wider initial wavepacket (Rq — 10). The results 
with A = 5 > A c 4.4 are shown in Fig. 5. One can find 
that the wavepacket is also trapped into several lattices 
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(at t « 100) as wavepacket moving. In this case, how- 
ever, the self-trapping occurs only during a certain time 
(t < 100, this time scale decreasing with A) but then the 
wavepacket start to diffuse. If we strengthen the nonlin- 
ear term to A = 20 >> A c = 4.4, the associated results 
are shown in Fig. 6. It indicates that the condensate is 
trapped almost in two optical lattices eventually and dif- 
fusion is not occur. Also, Fig. 7 provides the results with 
A < A c , as expected the diffusion result is obtained. 




FIG. 5: The dynamics of the 2D wavepacket with A = 5 > 
A c ~ 4.4, R = 10, po = 37r/4 



Figs. 8-9 show the numerical results for the dynamics of 
3D wavepacket. Differ from the results in 2D, when A > 
A c (the first panel in Fig. 8 and Fig. 9, respectively), the 
3D wavepacket gets into self-trapping state and finally, 
the wavepacket splits into two pieces and localized state 
appears. The results with A < A c (the second panel in 
Fig. 8 and Fig. 9, respectively) give the diffusion results. 
It is clearly that the above numerical results are in good 
agreement with our variational predications. It can be 
seen the dimension of the lattices has significant effect 
on the dynamic of BEC wavepackets moving in optical 
lattices. 



V. CONCLUSION 



To summarize, we have investigated the dynamics of 
a BEC wavepacket loaded into a deep multi-dimensional 




FIG. 6: The dynamics of the 2D wavepacket with A = 20 >> 
A c ~ 4.4, R = 10, po = 3tt/4 




FIG. 7: The dynamics of the 2D wavepacket with A = 3 < 
A c ~ 4.4, R = 10, po = 3tt/4 



periodic optical lattices both numerically and analyti- 
cally, focusing on the effect of the lattice dimension. Our 
study shows that the dynamics of the BEC wavepacket 
is quite different in ID, 2D or 3D systems. For exam- 
ple, the stable moving soliton and breather states that 
exist in ID lattice system no longer stand in 2D and 3D 
lattice systems. We also obtain an analytical criterion 
for the self-trapping state of BEC and phase diagram for 
self-trapping, soliton, breather, or diffusion of the BEC 
cloud for ID, 2D and 3D, respectively. The above ana- 
lytical results are confirmed by our directly solving the 
GPE. We hope that our studies will stimulate the exper- 
iments in the direction. 



t=0 




FIG. 8: The dynamics of the 3D wavepacket with A = 30 > 
A c ~ 18 (the first panel) and A = 15 < A c ~ 18 (the second 
panel), Ro = 2.5, po = 3n/4 
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t=20 



t = 200 




FIG. 9: The dynamics of the 3D wavepacket with A = 90 > 
A c ~ 79 (the first panel) and A = 50 < A c ~ 79 (the second 
panel), Ro = 10, po = 37r/4 
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